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ABSTRACT 

A Fourier technique for cloud motion estimation from pairs of pictures is 
described. This technique makes use of the phase of the cross spectral density; 
it allows motion estimates to be made for individual spatial frequencies, which 
are related to cloud pattern dimensions. Results obtained using this technique 
are presented and are compared with the results of a Fourier-domain cross- 
correlation scheme, using both artificial and real cloud data. It is concluded 
that this technique is relatively sensitive to the presence of mixtures of motions, 
changes in cloud shape, and edge effects. 
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A FOURIER APPROACH TO CLOUD MOTION ESTIMATION 
I. INTRODUCTION 

ml 

Sequences of high resolution photographs from satellites in geosynchronous 
orbit have been utilized by many investigators to measure cloud motion and de- 
rive winds. The techniques for extracting such information range from a 
manual analysis of closed "movie loops" (Fujita et al. , 1589; Hubert and 
Whitney, 1971) to completely automated methods in which the images and/or 
data are computer processed but selection of cloud targets is controlled manu- 
ally from an interactive display console (Serebreny et al. 1969; Smith and 
Phillips, 1972). 

The main advantage of manual and semi-automatic methods is in target 
selection, making it possible to incorporate a wide variety of criteria, not 
necessarily defined beforehand, for choosing cloud elements representative of 
atmospheric motion. These criteria are sometimes related to the size of the 
cloud element, the idea being that small elements move with the wind, while 
massive elements may disguise the actual motion of the atmosphere due to 
cloud formation and deformation processes. 

Early studies by Leese and Epstein (1963) have shown that one can identify 
the predominant dimensions of cloud patterns by a spectral analysis of satellite 
photographs. Subsequently, Leese et al. (1971) used the Fourier transform as 
an efficient computational tool for obtaining the cross-correlation between suc- 
cessive cloud cover pictures, from which cloud displacements are derived. 
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The fact that displacement information appears separately for each spectral 
component in the Fourier transform of the cross-co variance function was 
pointed out by Weinstein (1972). The combined advantages of the Fourier trans- 
form — as an efficient computational tool for measuring displacement and as a 
means for extracting cloud pattern dimensional information — have motivated 
the present study of the Fourier transform properties of real and simulated 
sequences of satellite cloud cover pictures and their cro c s -correlations. Spe- 
cifically, we describe an attempt to derive winds from cloud motion, while dis- 
criminating with respect to cloud size. 

The mathematical aspects of the Fourier transform method are outlined in 
Section 2; the results of a simulation study are presented in Section 3; a pre- 
liminary experiment with real pictures is described in Section 4; and finally, 
the outlook is discussed and conclusions drawn in Section 5. 

H. PROPERTIES OF THE FOURIER TRANSFORM 

A function defined in the space domain, say f(x,y), where x and y are co- 
ordinates in a two-dimensional Cartesian space, can be transformed into a 
function in the frequency domain by using the Fourier transform, F, as follows; 

F[f(x,y)] = f(x,y}e- Wta+v >»dxdy s F (u, v) (1) 


where i =~\/aT; u and v are frequency domain variables. Conversely, a function 
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in the frequency domain can be transformed into the space domain by using the 
inverse Fourier transform; 


F- 1 [F(u,v)j 



P(u,v) e 2wi ^ ,x+v ^ dudv a f(x,y) 


(2) 


For application to discrete data, die Fourier transform process is re- 
defined as 


F(u,v) 


f(x,y) 



N— 1 

-2tt[ 


l 

y=0 

f(x,y)e 

\M N j 

N— 1 

1 



2jri [ 

2j 

y=0 

F(u,v) e ' 

N/ 


(3) 


(4) 


where M and N are the sample sizes in the x and y directions, respectively. 

One indication of the relative displacement between successive pictures is 
by the location of the peak of the cross-covariance of the two images (Leese et 
ah , 1971). Letting f(x,y) and g(x,y) represent the brightness distribution of 
the two images, the cross-covariance is defined as 


C rg (x,y) 


-ff 


g(^-x,rj-y) d^clr? 


(5) 


The Fourier transform provides a fast way to calculate the cross-covariance 
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function. From equations (1) and (2), it follows that 

C f} ,(x,y) = F -1 [F !|! (u,v) - G(u,v)] (6) 

where F*(u, v) is the complex conjugate of F(u, v) while G(u, v) is the Fourier 
transform of g(x,y). The Fourier transform of C fu , namely F*G, is called Hie 
cross spectral density of f and g. The positions of the maxima of the cross- 
covariance function determine the relative displacement between the two pic- 
tures. The ratio of the computer times required for the Fourier transform 
method and Hie conventional method of calculating the cross -covariance is ap- 
proximately (N/2) N (Brigham and Morrow, 1967), where N is the linear size 
of the array. 

If the discrete version of the Fourier transform is used to compute the 
cross-covariance, as in (6), the resulting function is cyclic, which is equivalent 
to assuming that the functions in the space domain are periodic. In other words, 
the Fourier transform method assumes that when movement of objects takes 
place, what moves out of the picture at one boundary must enter at the opposite 
boundary. Further examination and discussion of this assumption is presented 
in Section m. 

The Fourier transform possesses still another property which makes it 
possible to take a different approach to the estimation of motion. Assume that 
uniform translation is the only difference between the two functions, i. e. , 
g(x,y) = f(x-a, y-b). Then from the definition of the Fourier transform in 
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equation (1) it follows (Briggs, 196 S) that 

G(u,v) = F(u,v) e“ 2ni( " ll+v b) (7) 

The Fourier transform at a specific frequency (u, v) can be represented by 
its amplitude and phase angle in the following manner: 

F(u,v) = lF(u,v) I e ip(u > v) (8) 

where |F| is the amplitude and p is the phase angle. Applying the same repre- 
sentation to G(u, v), we have 

G(u,v) = lF(u,v) I e'^'-^-^Kuii+vb) (9) 

Thus, 

F*(u,v) G(u.v) = lF(u,v)l 2 e- 2ff h«a+vb) (10) 

The term F* (u, v)G(u, v) is complex with amplitude |F! 2 and phase angle 
—2 7r (ua+vb), for all u's and v's. The differences between the phase angles of 
two neighboring frequencies have the following properties: 

ip(u,v) — ip(u + ] ,v) = 27ra (11) 

and 

sflLi.vl-tpfu.v+l) = 2?rb (12) 

We thus find that the difference in phase angle between neighboring terms in the 
F*G matrix yields the x and y components of displacement. 
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In the case where the two images differ only in that one is a uniform trans- 
lation of the other, the phase differences are the same between all neighboring 
frequencies. However, when the displacement between the two images is not 
uniform — and, specifically, when displacement is different for feature': with 
different characteristic dimensions — then the phase angle difference will be a 
function of frequency. In the following sections, we will examine to what extent 
the phase angle differences continue to yield estimates of the displacement and 
how these estimates are related to the size of features in the image. 

HI. SIMULATION STUDIES 

A number of simulation studies were conducted to investigate the compara- 
tive effectiveness of cross-correlation and phase difference techniques for es- 
timating motion {Lo and Parikh, 1973; Lo et al. , 1974). In particular, cases 
were studied that involved 

a. Distortions soich as rotation, scale change, brightness change, and 
noise, in addition to simple translation 

b. Mixtures of motion 

c. Edge effects: objects moving off one boundary without coming back 
in at the opposite boundary, so that the displacements are non-cyclic. 

It was found that if these conditions were not too strongly present, both tech- 
niques yielded reasonable estimates of the motion (or, in mixture cases, esti- 
mates of the predominant motion). However, as the distortion, mixture, or 
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edge factors become stronger, performance began to deteriorate. This situation 
is illustrated in the simulations shown in Figures 1 and 2. 

In Figure 1, a large mass of relatively dim cloud is moving at a compara- 
tively slow speed (displacement components X = 2 and Y = 1) relative to 
smaller clouds of higher brightness (X = 4, Y = 6). Figure la shows the two 
pictures, Figure lb shows their power spectra, and Figure 1c shows their 
cross-covariance, which has a peak corresponding to the displacement (4, 6) 
of the higher-contrast clouds. Figure Id shows the phase angles of the Fourier 
transform pi’oduct, and Figure le shows the resulting X and Y displacement 
estimates from the phase difference method (obtained b 3 <- taking horizontal and 
vertical differences in Figure Id). Figure If is a scatter plot of these esti- 
mates, which cluster around the displacement (4, 6). 

Figure 2 is analogous to Figure 1 for a pair of pictures (shown in Figure 2a) 
containing two sets of clouds with the same brightness but different velocities, 
represented by the displacements (6,0) and (2,2). Here the cross-covariance 
(Figure 2c) corresponds to the displacement (2,2) of the larger cloud mass; but 
the phase difference estimates (Figure 2e) show little tendency to cluster (Figure 2f). 

It would seem from these examples that flue phase difference estimation 
approach is more sensitive to the presence of mixtures of motions than the 
cross-covai’iance approach. Similar results can be obtained by studying simula- 
tions involving edge effects and other distortions. The details can be found in 
the two references cited at the beginning of the section. 
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IV. EXPERIMENTS WITH REAL DATA 

Several experiments with, pairs of real cloud windows were also carried 
out, as described by Lo et al. (1974). The windows used were selected from 
ATS-1 geosynchronous satellite images, taken at intervals of 47 minutes. Each 
window is 64 by 64 pixels, and has been scaled to the gray level range 0 to 63. 

The two windows shown in Figure 3a contain relatively small clouds. In 
order to better evaluate the synoptic environment in which these cloud patterns 
are embedded, two larger windows (256 by 256 pixels) containing the smaller 
windows at their centers are shown in Figure 3b. It is evident that the clouds 
in this scene tend to dissipate and re-form relatively fast, so that considerable 
distortion has occurred between the frames. Nevertheless, many of the cloud 
patterns in the first window are easily recognizable in the second. Based on 
such patterns, a hand estimate of the relative displacement was made using 
computer printouts of the pictures so that the displacement could be accu- 
rately measured. The estimated displacement, averaged over several cloud 
elements, Avas X = -5, Y = -3. 

The phase difference method applied to the pair of windows in Figure 3a 
gave rise to a displacement estimate of (-12, -4). This estimate was ob- 
tained as follows: X and Y displacement estimates were obtained for each 
spatial frequency component. The means and standard deviations of these es- 
timates were computed, and estimates deviating from the mean by more than 
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one standard deviation were discarded. The mean was then recomputed for the 
surviving estimates, and was used as the final estimate. 

The cross -covariance method, on the other hand, yielded a displacement es- 
timate of (-6, -3), as determined by the position of the cross-covariance peak. 

In tins example, therefore, the cross-covariance method yields better results 
{assuming the hand estimate to be most reliable). 

An attempt was made to improve the results obtained from the phase differ- 
ence method by applying it to selected bands of spatial frequencies, rather than 
to all frequencies. The following results were obtained: 


Band 

Frequency Range 

X 

Y 

a 

0 < yju 2 + v 2 < 5 

-8 

2 

b 

5 < V u 1 + v 2 <10 

-13 

16 

c 

10 < \u 2 + v 2 <64 

-12 

-4 


Note that the last result is the same as that obtained when all frequencies were 
used, and that those results still do not agree very well with the hand estimate. 
The cross-covariance estimates were also re-evaluated using only selected 
spatial frequencies {i. e. , the cross spectral density F*G was computed, the 
frequencies outside the selected band were suppressed, and the inverse Fourier 
transform of the result was computed to yield a "filtered" cross-covariance). 
Tlie results were as follows: 


T 


T 


T 


Frequency Band 

x X 

a 

-6 -3. 5 

b 

-5 -3 

c 

No clearcut peak 
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Results (a-b) are quite close to the result (-6, -3) obtained when no filtering was 
used, and also agree well with the hand estimate (-5, -3). 

V. DISCUSSION AND CONCLUSIONS 

The results of the experiments reported here indicate that both the cross- 
covariance and phase difference methods can be good estimators of motion when 
the objects being tracked do not change their shape, size and orientation to more 
than a limited degree. These teclmiques are less effective when a mixture of 
motions exists, unless one of the motions is strongly dominant. These prop- 
erties indicate that the Fourier transform phase difference estimation methods 
should be reliable in problems such as landmark matching, where the features 
do not change appreciably as viewed in satellite photographs. 

In the atmosphere, where the clouds grow, dissipate, rotate, and move 
relative to each other, the estimation of these motions using Fourier transform 
methods is not always reliable. 

The main obstacle to accurate cloud motion estimation is believed to be the 
presence of mixtures of motions. One possible way to solve this problem is to 
use spatial frequency filtering to separate the motions of different types of 
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cloud features, since cloud size is closely related on die one hand to spatial fre- 
quency and on the other hand to cloud motion. The effort to separate clouds by 
their sizes through frequency filtering is, however, thwarted by the fact that for 
real pictures each component in the frequency domain is a complicated function 
of the features in the space domain. Cloud patterns are not perfect sinusoidal 
functions. More than one frequency component is always needed to represent a 
cloud element, and conversely a high frequency component may receive contri- 
butions from small cloud elements , and also from sharp corners of a large 
cloud element. To further complicate the matter, the small cloud element or 
sharp corner may not exist in the second picture with which the first picture is 
being correlated. It is not found possible to effectively separate cloud types by 
spa dal frequency filtering. 

Another possible approach to separating' cloud types for motion estimation 
purposes is to segment die cloud cover windows using thresholding or pattern 
classification techniques (Lo and Mohr, 1974; Lo, 1975). The brightness, the 
equivalent black body temperature in various regions of the infrared spectrum, 
and the size of clouds are all closely related to the altitude of the clouds. Since 
clouds at certain altitudes in a limited region tend to move at the same velocity, 
thresholding techniques could be designed to separate clouds according to their 
brightness and/or spectral equivalent black body temperatures. Picture analy- 
sis techniques which determine cloud size should also be investigated for the 


separation of clouds. 


12 


It ca.n be concluded from this study that a technique wiiich separates clouds 
according to their motions must be designed and utilized to pre-process a pic- 
ture before applying Fourier transform motion estimation tecliniques such as 
the cross -covariance and phase difference methods. The cross-covariance 
method appears to be somewhat more reliable than the phase difference method, 
and should be preferred, even though the latter does present the advantages of 
providing multiple estimates and relative economy in computer time. 
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CAPTIONS 


Figure 1 

a. Gray levels of a pair of images simulating cloud motion; a large mass of 
relatively dim cloud undergoing a displacement (in pixels) of X = 2 and 

Y = 1, on which are superimposed smaller but brighter clouds under- 
going a displacement of X = 4, Y = 6. 

b. Power spectra of the images in a. The (0, 0) component is in the upper 
left corner, with wave number increasing downward and to the right. 

c. Cross-covariance of the images in a. Zero displacement corresponds to 
the upper left corner; X displacement increases to the right and Y dis- 
placement increases downward. 

d. Phase angles of the cross spectral density of the two images. The (0,0) 
component is in the upper left corner, with wave number increasing down- 
ward and to the right. 

e. X and Y displacement estimates derived from the cross spectral density. 

f. Scatter plot of the estimates in e. 

Figure 2 (a-f) 

Analogous to Figure 1, but for images simulating the motion of two sets of clouds 
with the same brightness but different velocities (i.e., displaced by different 
amounts) . 

Figure 3 

A pair of images selected as an example from the ATS-1 satellite, (a) The 
64 x 64 pixel portion of each image from which cloud motion is to be determined, 
(b) The surrounding area (256 by 256 pixels) in which the images in a are 
centrally imbedded. 
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Figure la. Gray levels of a pair of images simulating cloud motion: a large mass or relatively dim 
cloud undergoing a displacement (in pixels) of X = 2 and Y = 1, on which are superimposed smaller 
but brighter clouds undergoing a displacement of X — 4, Y = 6. 
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Figure lb. Power spectra of the images in a. The (0,0) component is in the upper 
left corner, with wave number increasing downward and to the right. 
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Figure le. X and Y displacement estimates derived from the cross spectral density. 
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Figure 2a. Analogous to Figure 1, but for images simulating the motion of two sets of clouds with 
the same brightness but different velocities (i.e., displaced by different amounts). 
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Figure 2b. Power spectra of the images in a. The (0,0) component is in the upper 
left corner, with wave number increasing downward and to the right. 
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Figure 2e. X and Y displacement estimates derived from the cross spectral density. 
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